package feap.compression_test_j2_std;

import java.util.ArrayList;

import math2.MathUtilities;
import fem2.AbstractStructuralStaticDemo;
import fem2.Constraint;
import fem2.Debugger;
import fem2.Element;
import fem2.Face;
import fem2.MaterialModel;
import fem2.Mesh;
import fem2.MeshGenerator;
import fem2.MeshPart;
import fem2.MeshUtilities;
import fem2.Model;
import fem2.Node;
import fem2.PatchManager;
import fem2.StringUtilities;
import fem2.element.StructuralElement;
import fem2.enu.State;
import fem2.enu.StressType;
import fem2.material.FluencyCriteria;
import fem2.material.LinearIsotropicHardeningRule;
import fem2.material.LinearKinematicHardeningRule;
import fem2.material.MLSrecoveryMM;
import fem2.material.SPRrecoveryMM;
import fem2.material.VonMisesMM;
import fem2.mesh_and_topology.IntegrationPointIndexer;
import fem2.mesh_and_topology.SpatialCriteria;
import fem2.mesh_and_topology.SphereSpatialCriteria;
import fem2.observer.GidStressObserver;
import fem2.observer.TimeStepObserver;
import fem2.pre_and_post.PatchMeshGenerator;
import fem2.strategies.Strategy;

public class Ivig02cc extends AbstractStructuralStaticDemo {

	PatchManager pm;
	IntegrationPointIndexer indexer;

	public Ivig02cc() {
		projectDir = "/home/hbui/kratos_janosch/fem2012-feap-validation/compression-test-j2-std";
		projectName = "Ivig02cc.gid";
		meshFileName = "mesh4x13";
	}

	@Override
	public Mesh createMesh() {
		int nx = 4;
		int ny = 13;
		double w = 0.04;
		double h = 0.14;
		MeshGenerator mg = new PatchMeshGenerator(nx + 1, ny + 1, 0, 0, w / nx, h / ny);
		return mg.getMesh(2);
	}

	@Override
	public Model createConditions(Model m) {
		Mesh mesh = m.getMesh();
		ArrayList<Node> bottom = MeshUtilities.seekNodesOnLine(mesh, 0, 0, 1, 0);
		ArrayList<Node> top = MeshUtilities.seekNodesOnLine(mesh, 0, 0.14, 1, 0);

		Constraint cx = new Constraint(false, true);
		MeshUtilities.seekNode(mesh, 0, 0).setConstraint(cx);

		Constraint cy = new Constraint(true, false);
		mesh.addConstraint(cy, bottom);

		Constraint cu = new Constraint(true, false);
		cu.setValue(1, -0.5e-5);
		mesh.addConstraint(cu, top);

		return m;
	}

	@Override
	public MaterialModel createMaterial(Model m) {
		double E = 0.2e5;
		double nu = 0.4;
		double t = 1.0;
		double rho = 1.0;
		// MaterialModel mm = new StVenantKirchhoffMM(E, nu, t, rho,
		// State.PLANE_STRAIN);

		double ft = 70;
		double K = -100;
		double H = 0;
		FluencyCriteria isoLaw = new LinearIsotropicHardeningRule(ft, K);
		FluencyCriteria kinLaw = new LinearKinematicHardeningRule(H);
		MaterialModel mm = new VonMisesMM(E, nu, t, rho, State.PLANE_STRAIN, isoLaw, kinLaw);

		// return mm;

		double r = 0.02;
		return new MLSrecoveryMM(indexer, new SPRrecoveryMM(pm, mm), r);
	}

	public MaterialModel createMaterial2(Model m) {
		double E = 0.2e5;
		double nu = 0.4;
		double t = 1.0;
		double rho = 1.0;
		// MaterialModel mm = new StVenantKirchhoffMM(E, nu, t, rho,
		// State.PLANE_STRAIN);

		double ft = 0.75 * 70;
		double K = -100;
		double H = 0;
		FluencyCriteria isoLaw = new LinearIsotropicHardeningRule(ft, K);
		FluencyCriteria kinLaw = new LinearKinematicHardeningRule(H);
		MaterialModel mm = new VonMisesMM(E, nu, t, rho, State.PLANE_STRAIN, isoLaw, kinLaw);

		// return mm;
		double r = 0.02;
		return new MLSrecoveryMM(indexer, new SPRrecoveryMM(pm, mm), r);
	}

	@Override
	public Element createElement(MeshPart mp, MaterialModel mm) {
		return new StructuralElement(mp, mm);
	}

	public Model createModel() {

		Mesh mesh = createMesh();

		Model m = new Model(mesh);

		pm = new PatchManager(mesh);
		indexer = new IntegrationPointIndexer(mesh);

		MaterialModel mm = createMaterial(m);
		MaterialModel mm2 = createMaterial2(m);

		SpatialCriteria sp = new SphereSpatialCriteria(0.035, 0.07, 0.001);

		for (int i = 0; i < mesh.countFaces(); i++) {
			Element e;
			Face f = mesh.getFace(i);
			if (sp.isSatisfied(f.computeCenter())) {
				e = createElement(f, mm2);
			} else {
				e = createElement(f, mm);
			}
			m.addElement(e);
		}

		m = createConditions(m);

		/*
		 * set the number of fields for structural problem
		 */
		mesh.setNumberOfFields(mesh.getDimension());

		StructuralElement dummy = (StructuralElement) m.getElement(0);
		for (int i = 0; i < dummy.getNumIntegrationPoints(); i++) {
			Debugger.watch(dummy.getIntegrationPoint(i));
		}

		return m;
	}

	@Override
	public void addObservers(Model m, Strategy s, TimeStepObserver o) {
		super.addObservers(m, s, o);

		if (meshFileName == null) {
			throw new Error("mesh file name was not set");
		}

		String meshName = StringUtilities.getFileName(meshFileName);

		s.addObserver(new GidStressObserver(m, o, projectDir, projectName, meshName,
				StressType.RECOVER_STRESS_ZZ));

		s.addObserver(new GidStressObserver(m, o, projectDir, projectName, meshName,
				StressType.RECOVER_STRESS_MLS));
	}

	public static void main(String[] args) {
		AbstractStructuralStaticDemo demo = new Ivig02cc();

		// double[] lambda = MathUtilities.createLambdaArray(1, 1, 2);
		double[] lambda = MathUtilities.createLambdaArray(1, 1, 420);

		demo.run(lambda);
	}
}
